Flowcytometer Biomass Calculation

Author

Rainer M Krug & Romana Limberger

Introduction

Details the calculation of the biomass and the problem occuring with FSC.A == 0

Code
library(LEEF.analysis)
library(LEEF.measurement.flowcytometer)
library(scattermore)

Calculate measurements

Code
if (!dir.exists(file.path(params$output_dir, "pre_processed"))){
  dir.create(file.path(params$output_dir, "pre_processed"), showWarnings = FALSE, recursive = TRUE)
  
  file.copy(
    file.path(params$pre_processed_dir, paste0("LEEF.fast.flowcytometer.", as.character(params$timestamp))),
    file.path(params$output_dir, "pre_processed"),
    recursive = TRUE,
    copy.mode = FALSE
  )
  
  file.rename(
    from = file.path(params$output, "pre_processed", paste0("LEEF.fast.flowcytometer.", as.character(params$timestamp))),
    to   = file.path(params$output, "pre_processed", "flowcytometer")
  )
}

if (!dir.exists(file.path(params$output_dir, "extracted_raw"))){
  dir.create(file.path(params$output_dir, "extracted_raw"), showWarnings = FALSE, recursive = TRUE)
  
  LEEF.measurement.flowcytometer::extractor_flowcytometer(
    input = file.path(params$output, "pre_processed"),
    output = file.path(params$output_dir, "extracted_raw"),
    raw = TRUE
  )
}

if (!dir.exists(file.path(params$output_dir, "extracted_log10"))){
  dir.create(file.path(params$output_dir, "extracted_log10"), showWarnings = FALSE, recursive = TRUE)
  
  LEEF.measurement.flowcytometer::extractor_flowcytometer(
    input = file.path(params$output, "pre_processed"),
    output = file.path(params$output_dir, "extracted_log10"),
    raw = FALSE
  )
}

Read Data

Size Standard

Code
ss <- data.frame(
  diameter_micrometer = c(1L, 3L, 10L), 
  mean_FSC.A = c(102045L, 641752L, 1959230L), 
  log_mean_FSC.A = c(5.008792, 5.807367, 6.292085)
)

Log10 transformed

Code
datadir <- file.path(params$output_dir, "extracted_log10", "flowcytometer")

# write.csv(
#   data.frame(
#     bacteria1 = c(0, 10000000, 10000000, 0),
#     bacteria2 = c(0, 0,        10000000, 10000000),
#     LNA = c(3,4.579784, NA, NA),
#     MNA = c(4.579785, 4.954243, NA, NA),
#     HNA = c(4.954244, 7, NA, NA),
#     algae1 = c(3, 3, 7, 7),
#     algae2 = c(4, 7, 7, 4)
#   ),
#   file = file.path(datadir, "gates_coordinates.csv")
# )

traits_log10 <- LEEF.measurement.flowcytometer::extract_traits(
  input = datadir,
  particles = params$particles,
  metadata_flowcytometer = read.csv(file.path(datadir, "metadata_flowcytometer.csv"))
)
{"timestamp": "2023-03-06T15:48:42+0100", "log_lvl": "INFO", "log_msg": "   reading data  ..."}
   reading data  ...
{"timestamp": "2023-03-06T15:48:57+0100", "log_lvl": "INFO", "log_msg": "   done"}
   done
{"timestamp": "2023-03-06T15:48:57+0100", "log_lvl": "INFO", "log_msg": "########################################################"}
########################################################
Code
paste0("Count       : ",  traits_log10[[1]]$FSC.A |> length())
[1] "Count       : 5048859"
Code
paste0("Range       : ",  paste0(traits_log10[[1]]$FSC.A |> range(), collapse = " -- "))
[1] "Range       : 0 -- 7.22471987004958"
Code
paste0("Portion == 0: ", ((traits_log10[[1]]$FSC.A == 0) |> sum()) / (traits_log10[[1]]$FSC.A |> length()))
[1] "Portion == 0: 0.127955444982718"
Code
traits_log10[[1]]$FSC.A |> table() |> head()

                0 0.301029995663981 0.477121254719662 0.602059991327962 
           646029              1456              1440              2799 
0.698970004336019 0.778151250383644 
             1464              1421 

Raw

Code
datadir <- file.path(params$output_dir, "extracted_raw", "flowcytometer")

# write.csv(
#   x = data.frame(
#     bacteria1 = c(0, 15000000, 15000000, 0   ),
#     bacteria2 = c(0, 0,    15000000, 15000000),
#     LNA = c(3,4.579784, NA, NA),
#     MNA = c(4.579785, 4.954243, NA, NA),
#     HNA = c(4.954244, 7, NA, NA),
#     algae1 = c(3, 3, 7, 7),
#     algae2 = c(4, 7, 7, 4)
#   ),
#   file = file.path(datadir, "gates_coordinates.csv")
# )


traits_raw <- LEEF.measurement.flowcytometer::extract_traits(
  input = datadir,
  particles = params$particles,
  metadata_flowcytometer = read.csv(file.path(datadir, "metadata_flowcytometer.csv"))
)
{"timestamp": "2023-03-06T15:49:00+0100", "log_lvl": "INFO", "log_msg": "   reading data  ..."}
   reading data  ...
{"timestamp": "2023-03-06T15:49:13+0100", "log_lvl": "INFO", "log_msg": "   done"}
   done
{"timestamp": "2023-03-06T15:49:13+0100", "log_lvl": "INFO", "log_msg": "########################################################"}
########################################################
Code
paste0("Count       : ",  traits_raw[[1]]$FSC.A |> length())
[1] "Count       : 5048859"
Code
paste0("Range       : ",  paste0(traits_raw[[1]]$FSC.A |> range(), collapse = " -- "))
[1] "Range       : 0 -- 16777215"
Code
paste0("Portion == 0: ", ((traits_raw[[1]]$FSC.A == 0) |> sum()) / (traits_raw[[1]]$FSC.A |> length()))
[1] "Portion == 0: 0.127662903638228"
Code
traits_raw[[1]]$FSC.A |> table() |> head()

     0      1      2      3      4      5 
644552   1477   1456   1440   2799   1464 

Some Plots

Code
i0 <- traits_raw[[1]]$FSC.A == 0
i1 <- traits_raw[[1]]$FSC.A == 1

The Size standard and the extrapolation

Code
plot(ss$mean_FSC.A, ss$diameter_micrometer, xlim = c(0, 3000000), ylim = c(0, 20))
abline(a = params$length_intercept, b = params$length_slope)
abline(v = min(traits_raw[[1]]$FSC.A))
abline(v = max(traits_raw[[1]]$FSC.A))
abline(h = 0)

Scatterplots

FSC-A auf der x-Achse und SSC-A auf der y-Achse

Code
scattermoreplot(traits_raw[[1]]$FSC.A, traits_raw[[1]]$SSC.A, xlab = "FSC.A", ylab = "SSC.A")

Code
x <- traits_log10[[1]]$FSC.H
x[x==0] <- 1
x <- log10(x)
scattermoreplot(traits_log10[[1]]$FSC.A, x, xlab = "log10(FSC.A)", ylab = "log10(SSC.A)")

FSC-A auf der x-Achse und FL1-A auf der y-Achse

Code
scattermoreplot(traits_raw[[1]]$FSC.A, traits_raw[[1]]$FL1.A, xlab = "FSC.A", ylab = "FL1.A")

Code
x <- traits_log10[[1]]$FL1.A
x[x==0] <- 1
x <- log10(x)
scattermoreplot(traits_log10[[1]]$FSC.A, x, xlab = "log10(FSC.A)", ylab = "log10(FL1.A)")

FSC-A auf der x-Achse und FL3-A auf der y-Achse

Code
scattermoreplot(traits_raw[[1]]$FSC.A, traits_raw[[1]]$FL3.A, xlab = "FSC.A", ylab = "FL3.A")

Code
x <- traits_log10[[1]]$FL3.A
x[x==0] <- 1
x <- log10(x)
scattermoreplot(traits_log10[[1]]$FSC.A, x, xlab = "log10(FSC.A)", ylab = "log10(FL3.A)")

FSC-A auf der x-Achse und FSC-H auf der y-Achse

Code
scattermoreplot(traits_raw[[1]]$FSC.A, traits_raw[[1]]$FSC.H, xlab = "FSC.A", ylab = "FSC.H")

Code
x <- traits_log10[[1]]$FSC.H
x[x==0] <- 1
x <- log10(x)
scattermoreplot(traits_log10[[1]]$FSC.A, x, xlab = "log10(FSC.A)", ylab = "log10(FSC.H)")

And some more density plots

Raw FSC.A and FSC.H

Code
traits_raw[[1]]$FSC.A |> density(bw = 10, to = 5000) |> plot(main = "FSC.A")

Code
traits_raw[[1]]$FSC.H |> density(bw = 10, to = 10000) |> plot(main = "FSC.H")

log10 FSC.A and FSC.H

Code
traits_log10[[1]]$FSC.A |> density(bw = 0.02) |> plot(main = "log10(FSC.A)")

Code
x <- traits_log10[[1]]$FSC.H
x[x==0] <- 1
x |> log10() |> density(bw = 0.02) |> plot(main = "log10(FSC.H)")

  • BLACK Line: All particles with FSC.A == 0 (will tbe transformed to 0 using log10)
  • RED Line: All particles with FSC.A == 1 (will tbe transformed to 0 using log10)
  • Black Line: All particles with FSC.A != 0

SSC.A

Code
traits_raw[[1]]$SSC.A[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$SSC.A[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$SSC.A[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL1.A

Code
traits_raw[[1]]$FL1.A[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL1.A[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL1.A[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL2.A

Code
traits_raw[[1]]$FL2.A[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL2.A[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL2.A[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL3.A

Code
traits_raw[[1]]$FL3.A[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL3.A[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL3.A[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL4.A

Code
traits_raw[[1]]$FL4.A[i0] |> density(to = 2000, bw = 1) |> plot()
traits_raw[[1]]$FL4.A[i1] |> density(to = 2000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL4.A[!i0] |> density(to = 2000, bw = 1) |> lines(col = "green")

FL2.A

Code
traits_raw[[1]]$FL2.A[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL2.A[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL2.A[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FSC.H

Code
traits_raw[[1]]$FSC.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FSC.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FSC.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

SSC.H

Code
traits_raw[[1]]$SSC.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$SSC.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$SSC.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL1.H

Code
traits_raw[[1]]$FL1.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL1.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL1.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL2.H

Code
traits_raw[[1]]$FL2.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL2.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL2.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL3.H

Code
traits_raw[[1]]$FL3.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL3.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL3.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

FL4.H

Code
traits_raw[[1]]$FL4.H[i0] |> density(to = 10000, bw = 1) |> plot()
traits_raw[[1]]$FL4.H[i1] |> density(to = 10000, bw = 1) |> lines(col = "red")
traits_raw[[1]]$FL4.H[!i0] |> density(to = 10000, bw = 1) |> lines(col = "green")

Width

Code
traits_raw[[1]]$Width[i0] |> density(to = 500, bw = 1) |> plot()
traits_raw[[1]]$Width[i1] |> density(to = 500, bw = 1) |> lines(col = "red")
traits_raw[[1]]$Width[!i0] |> density(to = 500, bw = 1) |> lines(col = "green")

Now to the sizes

Code
#|

i <- traits_log10[[1]]$FSC.A != 0

traits_log10[[1]]$length[i] <- 10^(traits_log10[[1]]$FSC.A[i]) * params$length_slope + params$length_intercept
traits_log10[[1]]$volume <- 4/9 * pi * traits_log10[[1]]$length^3

traits_log10[[1]]$length |> density(na.rm = TRUE) |> plot(main = "length in micro meter")

Code
traits_log10[[1]]$length |> density(na.rm = TRUE, from = 1) |> plot(main = "length in micro meter (l > 1")

Code
traits_log10[[1]]$length |> density(na.rm = TRUE, from = 1, to = 30) |> plot(main = "length in micro meter (1 < l < 30_")